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ABSTRACT 



We present the results of a series of axisymmetric time- dependent 
magnetohydrodynamic (MHD) simulations of the propagation of cooling, 
overdense jets. Our numerical models are motivated by the properties of 
outflows associated with young stellar objects. A variety of initial field strengths 
and configurations are explored for both steady and time-variable (pulsed) jets. 
For the parameters of protostellar jets adopted here, even apparently weak 
magnetic fields with strengths B > 60/zG in the pre-shocked jet beam can 
have a significant effect on the dynamics, for example by altering the density, 
width, and fragmentation of thin shells formed by cooling gas. Strong toroidal 
fields (> 100/iG) with a radial profile that peaks near the surface of the jet 
result in the accumulation of dense shocked gas in a "nose cone" at the head 
of jet. We suggest that this structure is unstable in three-dimensions. A linear 
analysis predicts that axisymmetric pinch modes of the MHD Kelvin-Helmholtz 
instability should grow only slowly for the highly supermagnetosonic jets studied 
here; we find no evidence for them in our simulations. Some of our models 
appear unstable to current-driven pinch modes, however the resulting pressure 
and density variations induced in the jet beam are not large, making this 
mechanism an unlikely source of emission knots in the jet beam. In the case 
of pulsed jets, radial hoop stresses confine shocked jet material in the pulses 
to the axis, resulting in a higher density in the pulses in comparison to purely 
hydrodynamic models. In addition, if the magnetic field strength varies with 
radius, significant radial structure is produced in the pulses (the density is 
strongly axially peaked, for example) even if the density and velocity in the jet 
follow a constant "top-hat" profile initially. 



Subject headings: hydrodynamics — magnetohydrodynamics — ISM: jets and 
outflows — galaxies: jets 
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1. Introduction 

The most promising mechanism for the production of supersonic, highly collimated jets 
from low mass young stellar objects is by magnetic forces associated either directly with 
an accretion disk (Konigl & Pudritz 1999), or with the interaction between an accretion 
disk and a magnetized central star (Shu et al. 1999). Magnetic fields are also thought to 
contribute to the collimation of the jets on larger scales, (although only slowly so that the 
observed jet may only be a part of a much broader wind, e.g., Ostriker 1997). Thus, unless 
the outflowing material is highly resistive (which seems unlikely), protostellar jets should 
contain a dynamically important magnetic field which may affect both the propagation 
and stability of the outflow. Observation of the magnetic field strength associated with 
protostellar jets is difficult. However, by fitting one-dimensional radiative shock models to 
the observed line ratios in the bow shock of HH47, Morse et al. (1993) inferred an upper 
limit to the magnetic field in the ambient medium upstream of the jet of ~ 30/iG, a value 
which they argued was too small to affect the dynamics except by increasing the cooling 
lengths behind radiative shocks. More recently, Ray et al. (1997) have reported direct 
evidence for strong fields in the outflow associated with the source T Tau S through the 
detection of opposite circular polarization in the two spatially resolved outflow lobes. These 
authors infer very high field strengths within the lobes; several Gauss at a distance of few 
tens of AU from the central star. 

If protostellar jets contain strong magnetic fields as expected from theory, there may 
be a signature of such fields in their dynamics. While hydrodynamical studies of cooling, 
dense protostellar jets are widely available in the literature (e.g., see Raga 1995; Cabrit 
1997 for reviews), MHD models are less well explored. MHD studies of extragalactic (i.e., 
underdense and non-radiative) jets have been reported by Clarke, Norman, & Burns (1986), 
and Lind et al. (1989) for the case of toroidal fields, and by Kossl et al. (1990a; 1990b) for 
the case of toroidal and axial fields (see also Clarke 1996 for a recent review). For strong 
fields, characterized by a small value of (3 = 8irP/B 2 where B is the field strength and 
P the thermal pressure (although note such fields may still have an energy density small 
compared to the kinetic energy density in the flow), it is found that the cocoon formed by 
lateral expansion of hot, shocked gas from the head of the jet is strongly inhibited. Instead, 
a magnetically confined "nose cone" of shock processed jet material is formed between 
the bow shock and Mach disk. In addition, the stability properties of the jet beam itself 
are strongly altered by the presence of a magnetic field. MHD studies of the propagation 
of overdense, non-radiative jets into a uniform ambient medium which contains a helical 
magnetic field everywhere have been reported by Todo et al. (1992). They describe results 
for a number of field strengths and approximate the effects of radiative cooling with an 
equation of state with 7 = 1.2. The formation of a nose cone, as well as suppression of the 
cocoon, is evident for strong fields. In subsequent three-dimensional simulations (Todo et 
al. 1993), it was found that a magnetic field strength of 70/zG was sufficient to produce 
nonaxisymmetric kink instabilities in the jet beam which the authors suggest may be 
relevant to bends observed in the jets from some protostellar systems. This field strength is 
approaching the value inferred by Morse et al. (1993). 

More recently, several authors have reported results from MHD studies of dense jets 
in which optically thin radiative cooling is directly coupled to the dynamics. Frank et 
al. (1998, hereafter FRJN-C) report axisymmetric studies which confirm the formation of 
nose cones, and (for the particular form of the initial field distribution adopted in the jet) 
the presence of pinch modes driven by magnetic tension. Cerqueira et al. (1997, hereafter 
CGH; see also de Gouveia Dal Pino & Cerqueira 1996) have used an SPH code which 
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has been extended with an algorithm designed to represent magnetic forces to study the 
propagation of dense, cooling jets in three-dimensions. They find that the fragmentation of 
the dense shell formed by cooling at the head of the jet is strongly affected by magnetic 
forces, and they use these results to conclude that the clumpy structure observed in jets 
provides evidence that near the head of the jet the field must be primarily axial rather than 
helical. In their latest simulations (Cerqueira & de Gouveia Dal Pino 1999, hereafter CG), 
these authors also report evidence for both modes of the MHD Kelvin-Helmholtz (KH) 
and a current-driven pinch instability in models of protostellar jets. Finally, Frank et al. 
(1999a; see also Lery & Frank 1999) have studied the propagation of MHD jets launched by 
Keplerian accretion disks using analytic models for magneto-centrifugal outflows to describe 
the jet structure at the nozzle. A variety of complex behavior is seen in this case, some 
of which results from MHD effects and some which results from the non-uniform radial 
structure of the jet beam. 

Given the apparent sensitivity of jet dynamics to the assumed initial orientation and 
strength of the magnetic field within the jet beam suggested by published numerical results, 
it is clearly useful to survey different initial field profiles. In this paper, we present the 
results from a large number of two-dimensional, axisymmetric MHD simulations of dense, 
cooling jets. To emphasize the effects produced solely by the magnetic field, we study jets 
with a uniform ("top-hat") radial profile propagating into an initially uniform ambient 
medium. Since the theory of the production of such jets suggests a primarily toroidal 
field within the outflow far from the source, we focus on toroidal fields with different 
strengths and radial profiles. In most of our models, the ambient medium into which the 
jet propagates is unmagnetized, although in some cases we add a poloidal field which 
threads both the jet and ambient gas initially. Given the importance of nonaxisymmetric 
instabilities to the dynamics of magnetized jets (Todo et al. 1993), fully three-dimensional 
MHD simulations of cooling jets are important (CGH; CG), however in order to cover 
parameter space efficiently we confine the models described in this study to axisymmetry. 
Moreover, for the highly supermagnetosonic jets studied here, non-axisymmetric modes 
of, e.g., the K-H instability grow slowly and will not affect the jet in regions close to the 
nozzle. Of course, precession of the jet beam, a clumpy ambient medium, or a much lower 
magnetosonic Mach number could introduce three-dimensional effects much earlier. Our 
results extend previous studies in that we evolve the jet for much longer (using a grid 
roughly three times larger than used by FRJN-C or CG), and we focus our attention on 
strongly magnetized jet beams propagating into an unmagnetized ambient medium. Where 
there is overlap, we compare and contrast our results to previous studies throughout this 
paper. 

The paper is organized as follows. We describe our numerical techniques in §2. As an 
aid to the interpretation of our results, we review the relevant stability properties of MHD 
jets in §3. Our results for steady are presented in §4, and for pulsed jets in §5. In §6 we 
summarize and conclude. 



2. Method 
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2.1. MHD Algorithms 

The dynamical evolution of a non-relativistic magnetized jet is given by solutions to 
the equations of MHD: 

§ + pV.v = 0, (1) 

Dv 

p— = - Vp + (Vx B x B) /4tt (2) 

fJL*=-pV.<v-A + H, (3) 
iJr p 

^ = Vx(vxB). (4) 

In these equations, D/Dt = d/dt + v • V is the convective derivative, A and H represent 
cooling and heating terms, respectively, and the other symbols have their usual meanings. 
We assume an ideal gas law for the equation of state; that is, p = (7 — l)e, where 7 = 5/3, 
and e is the internal energy density. We adopt the assumptions of ideal MHD in writing 
equations (1) - (4), in particular we assume the magnetic field is perfectly coupled to the 
fluid. While temperatures in the jet beam and in the strong shocks at the head of the jet 
are probably high enough (and therefore the gas is ionized enough) to keep the field well 
coupled on sub-parsec scales (Frank et al 1999b), it is possible that in slower shocks near 
the wings of the outflow, non-ideal MHD effects may become important, especially in dense 
ambient gas. However, the study of non-ideal MHD effects (for example, the structure of 
"C-type" shocks, Draine & McKee 1993) in the outflow is beyond the scope of this paper. 

The second and third terms on the right-hand-side of equation (3) represent energy 
loss via optically thin cooling and heating respectively. An accurate representation of these 
terms is one of the most important, and yet challenging, aspects of dynamical studies of 
cooling jets. Previous work has shown significant differences in the evolution of cooling 
jets using coronal versus non-equilibrium cooling rates (Stone & Norman 1993). Recently, 
Raga et al. (1997) have developed an accurate and efficient formulation to approximate 
non-equilibrium cooling rates in interstellar gas. However, given that the focus of this 
paper is the effect of magnetic fields on the propagation of cooling jets, we have chosen 
to adopt only a simple formulation for the optically thin cooling rate: the coronal cooling 
curve for interstellar gas calculated by Dalgarno & McCray (1972). We allow gas to cool to 
a temperature of 100 K, below this we set the cooling rate to zero. 

We use the ZEUS code (Stone & Norman 1992a; 1992b) to solve the MHD equations 
in cylindrical coordinates assuming axisymmetry. The stiff cooling and heating terms are 
differenced implicitly, resulting in a nonlinear equation at every grid point which is solved 
using a Newton-Raphson iteration scheme. This step is operator split from the rest of the 
MHD equations. No other additions or modifications to the ZEUS algorithm were required. 



2.2. Initial and Boundary Conditions 

The axisymmetric simulations are performed on a grid of size < r < 20Rj in the 
transverse direction, and < z < 100-Rj in the axial direction, where Rj is the initial 



- 6- 



jet radius. Except where indicated we use 400 uniform zones in the transverse direction 
and 2000 uniform zones in the axial direction, so that 20 zones span the initial jet radius. 
Reflecting boundary conditions are used along the inner r boundary, while outflow boundary 
conditions are used along the outer r, outer z, and inner z boundary for r > Rj. For r < Rj 
along the inner z boundary, inflow boundary conditions are used with the variables held 
fixed at the values appropriate to the initial equilibrium structure of the jet. 

Initially the jet is assumed to be perfectly collimated. Motivated by observations of 
various protostellar jets (e.g., Ray 1996), we adopt the following values as characteristic of 
observed systems: an initial jet radius Rj = 2.5 x 10 15 cm, a uniform density rij = 1000 cm -3 , 
and a uniform jet velocity Vj = 332 km s _1 . The gas pressure at the surface of the jet is 
taken to be po = 1.38 x 10~ 10 dyne cm -2 , corresponding to a temperature of Tj = 10 3 K 
and sound speed aj = 5.25 km s -1 . As described below, in order to ensure magnetostatic 
equilibrium in the radial direction, in some simulations we adopt a radially varying pressure 
profile. The ambient medium is stationary, with uniform density n a = 100 cm -3 (so that 
the jet is overdense, i.e. rj = rij/n a = 10), and temperature T a = 50 K. The sound speed in 
the ambient gas is a a — 1.17 km s -1 . Note that with these choices, the jet is overpressured 
with respect to the ambient medium, i.e., p a = po/a and a = 200. Several purely adiabatic 
simulations were performed for comparison with the cooling models, in this case we have 
taken p a = po, T a = 10,000 K and a a = 16.6 km s -1 , and the jet is initially in pressure 
equilibrium with the ambient medium. Physically our radiatively cooling models represent 
a dense jet beam surrounded by a hotter cocoon embedded in a cold ISM. The kinetic model 
of a dense jet embedded in a less dense cocoon which is in turn ensheathed by a dense and 
cold molecular gas is consistent with detailed observation of a number of protostellar jet 
systems including HH111 (Nagar et al. 1997). 

Since the strength and orientation of the magnetic field in protostellar jets is not 
determined observationally, there is considerable freedom in specifying the initial field 
configuration. We choose to focus on primarily toroidal fields in this paper (although 
we do present models which include both toroidal and axial fields) as theoretical models 
of magnetically driven outflows predict the field should become primarily toroidal 
asymptotically. In the absence of a definitive theory for the radial profile of the magnetic 
field, we assume the toroidal magnetic field is zero on the jet axis, achieves a maximum at 
some radial position in the jet interior r m , and returns to zero at the outer jet boundary. 
This implies that all currents which support the field are confined within the jet (although 
as the simulation progresses magnetic field, and therefore current, is convected into the 
cocoon). Within the jet, we adopt the simple profile for the toroidal magnetic field 



BJr) 



B^mf- < r < r m 

Rj-r 

'i m Rj-r m ' m — ' — 



B^m^Eh r m <r<R 1 (5) 



Rj < r 



This profile is identical to that used by Lind et al. (1989) and by FRJN-C when 
< r < r m and when r > Rj but differs for r m < r < Rj where they used a force-free 
B<p(r) = B^ m {r m /r) (which requires a large return current located at the jet surface). CG 
have studied both longitudinal and helical fields of the form used by Todo et al. (1992) 
throughout the computational domain. Our toroidal magnetic field corresponds to a uniform 
current density inside radius r m , and a distributed return current in the outer portion of the 
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jet. In adiabatic and cooling jet simulations designed to be compared directly, the toroidal 
magnetic field reaches a maximum value at r m /Rj = 0.9 with /3i = f3 rn = 87rpQ/B^ m = 1 
(Bfam = 58.8 jiG) and with /9i/ 4 = f3 m — 0.25 {B^ m = 117.6 fiG). To illustrate the toroidal 
magnetic field profile adopted here, we plot equation (5) for j3 m = 1 and r m = 0.9 in Figure 
1. 



The equation of hydromagnetic equilibrium 



Pj(r) + 




4irr 



(6) 



has been used to establish a suitable radial gas pressure profile in the jet 
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(7) 



(8) 



The radial pressure profile which gives a magnetostatic equilibrium is shown in Figure 1 
where we plot equation (7) for j3 m = 1 and r m = 0.9. Since the jet density is held constant, 
Tj has the same functional form as pj. In some of our simulations with strong toroidal fields, 
the condition of magnetohydrostatic equilibrium within the jet beam requires negative 
pressures in the neighborhood of r m . In these cases we simply fix p m = 0, and scale the 
pressure profile given by equation (7) by an appropriate amount so that this condition is 
met. Of course, this means that the jet beam is no longer in strict radial equilibrium. 

The profile adopted for the magnetic field is an important difference between this 
work and previous studies (FRJN-C, CG). To allow us to study the effect of the field 
geometry, we have also performed a few simulations in which the field contains a uniform 
axial component as well as the toroidal distribution given above, making the field helical. 
We will also describe the evolution of force-free fields with B^ oc 1/r. We find that the 
detailed structure of a magnetized radiative jet is sensitive to the precise field geometry, 
thus uncertainty about this geometry is a serious limitation to modeling observed systems. 

To construct models of pulsed jets, we apply a purely sinusoidal perturbation in time 
so the axial component of the velocity is 



u z {x = 0) = Vj(l + Asinut) , 



where we choose A 



(9) 

16.6 km s _1 



0.25, and for all of the simulations oj = 25a /Rj where a 
is a fiducial reference speed (equivalent to the sound speed at a temperature of 10 4 K). 

The important dimensionless parameters which characterize the evolution of the jet 
are the density ratio rj (taken to be 10 for all models presented here), and the internal 
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magnetosonic Mach number M^ ms = Vj/a rns , where a rns = (a 2 + V^) 1 / 2 and Va is the Alfven 
speed. In simulations without toroidal magnetic field or with force-free toroidal magnetic 
field the pressure, density, and temperature are uniform across the jet, so that the internal 
sonic Mach number Mj = Vj/aj is constant. However, a radially varying magnetic pressure 
requires a radially varying thermal pressure in order to initialize a jet beam in radial force 
balance and the sonic and magnetosonic Mach numbers are not constant. In general, we 
characterize all our models by the magnetosonic Mach number on the jet axis. However, 
since internal dynamics and timescales involve wave propagation across the jet through a 
medium of varying magnetosonic speed we also define linear radial averages as, for example, 

Mj,ms = 4" I*' M hms {r)dr. (10) 
rij Jo 



2.3. Tests of the Method 

Tests of the MHD algorithms used here have been reported in a number of papers 
(Stone & Norman 1992a, 1992b; Stone et al. 1992; Hawley & Stone 1995). We have also 
confirmed that we are able to reproduce the results of previous MHD studies of jets. In 
particular, we have performed simulations of adiabatic, underdense, magnetized jets using 
parameters identical to Clarke, Norman, & Burns (1986); we find excellent agreement 
between our results with those reported by these authors. In particular, we find production 
of a nose cone with a similar aspect ratio, strong pinch modes within the jet beam due to 
lack of initial radial equilibrium, and strong suppression of vortical motion in the cocoon by 
magnetic stresses. 



3. Jet stability 

Before describing the results of our simulations, let us review the basic stability 
properties of supermagnetosonic jets in order to aid interpretation. The presence of toroidal 
magnetic fields and the poloidal current can give rise to current driven pinch modes in 
addition to velocity shear driven Kelvin-Helmholtz pinch modes. In the absence of any 
zero order velocity shear, e.g., well inside the velocity shear surface at the jet boundary, 
Begelman (1998) has shown that a beam containing a toroidal magnetic field can be 
unstable to current driven pinch modes when 

dlnB T/3(r)-2 

cilnr T/3(r) + 2 1 1 

where the plasma beta f3(r) = 8np(r) / B 2 (r) , and T is the adiabatic index. In general, 
shorter pinch lengths will grow faster. If the beam contains a toroidal magnetic field and a 
weaker uniform axial magnetic field the condition becomes 

dlnB 1. , B z n2 

-— > -1 + - (kr—f- 2 , 12 
Gtlnr 2 B^r) 
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where k is the wavenumber of the mode; in general, shorter pinch lengths, i.e., larger k, are 
stabilized. 

The velocity driven Kelvin-Helmholtz unstable pinch body modes might appear near 
resonances which, in the absence of radial variation in parameters in the jet and in the 
external medium, are given by (e.g., Hardee, Clarke, & Rosen 1997; Hardee & Stone 1997) 

uRj/a c w (2m + 1/2)tt/2 (12a) 

X/Rj » -J^ (12b) 

where u, m, and A are the angular frequency, body mode number m — 1, 2, 3, etc., and 
wavelength, respectively. In the expressions above M c = t>j/a c and a c is the sound speed in 
the cocoon medium immediately outside the jet and not the sound speed in the undisturbed 
ambient medium. At these resonances the pinch modes grow with spatial growth rate 

k T w -(2M jjms R j )~ 1 In {4wRj/a c ) , (13) 



where kj is the imaginary part of the wavenumber. The spatial growth length over which 

a linear perturbation increases in amplitude by a factor e is £ = kj 1 ~ Mj tTns Rj. On the 

high Mach number axisymmetric flows that we simulate here these Kelvin-Helmholtz modes 
grow relatively slowly. 



4. Steady Jet Simulations 

Simulation parameters for the steady jet models discussed here are listed in Table 
1. In what follows, we report times in units of a sound crossing time of the jet beam, 
t = Rj/dj rs 152 yr where dj = 5.25 km s _1 . Magnetosonic crossing times of the jet 
beam are in the range, r ms = Rj/a ms pa 150 — 58 yr where a ms is the linear average of 
the magnetosonic speed. A grid crossing time can be defined by r g = L gr id/vh ~ 314 yr, 

where Vh = vj ^ 250 km s 1 is the velocity of advance of the head of a 

cold (thermal and magnetic pressure much less than the ram pressure) jet. We stop the 
evolution when the bow shock has reached ~ 100-Rj. The spatial evolution of the jet's 
transverse structure is then evaluated at an axial distance of 68.5Rj (a location behind 
the region strongly affected by jet head vortices), and the spatial evolution of the cocoon's 
structure can be evaluated through comparison of cocoon properties at this distance and 
near to the inlet. 



4.1. Magnetized Adiabatic Jets 

In order to understand the dynamics of our magnetized radiative jets, it is necessary 
to first compute non-radiative jet models with an identical magnetic field profile. We 
have computed three models: a purely hydrodynamical jet ((3^?), a weak field jet ((3^ d ), 
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and a strong field jet (Pi/ d 4 ). Table 1 gives additional parameter values adopted for each 

simulation. Unlike our radiatively cooling jet models, all three adiabatic models are initially 
in pressure equilibrium with the ambient medium; i.e. we adopt T a = 10 4 K rather than 
the T a = 50 K used for the cooling jets. Figure 2 plots the logarithm of the density at the 
end of each of these three simulations. Figure 3 plots the toroidal magnetic field strength 
at the same times for the two magnetized jet models. Note only the final 50Rj of the 
computational domain is shown in both figures for clarity. 

Comparison of the final structure of each jet shows remarkably little difference between 
these models. In the purely hydrodynamical jet simulation P^f the jet remains well 
collimated within its initial radius inside a low density cocoon of typical density one-half 
the ambient and one-twentieth the jet density, i.e., n c 50 cm -3 , and has a typical radius 
« 10.R,-. 

The weak toroidally magnetized jet simulation, @f d , experiences some slight collimation 
from the toroidal magnetic field so that by an axial distance of 68.5 Rj, where ij ^ Rj, 
the toroidal field maximum has increased by about 7% as it has moved inwards with 
the contracted jet. However, unlike the P^f jet we find three significant pressure pulses 
with amplitudes up to twice Pj(r) and spacing of about 15Rj (with the first pulse at an 
axial distance of ~ 30Rj) in the jet interior when 0.1 < r/rj < 0.6. The amplitudes are 
largest at r ss 0.3rj. This location is not consistent with pressure perturbations produced 
by the Kelvin-Helmholtz instability at the jet surface. The confinement of these pressure 
perturbations to the jet interior is consistent with a current driven instability, moreover this 
region satisfies the condition for instability (eq. [11]) found by Begelman (1998). 

The strong toroidally magnetized jet simulation, Pyj^, cannot be introduced at the inlet 

with a true equilibrium pressure distribution: this would require a negative minimum jet 
pressure. As a result the jet spine expands to rj m 1.5Rj at an axial distance of Q8.5Rj and 
then recollimates somewhat under the influence of higher external pressure near to the jet 
head and the stronger toroidal magnetic field. At axial distances between 80 — 95Rj there is 
evidence for pressure pulses associated with pinching in the jet interior at r < 0.2rj. Three 
main pressure pulses with spacing of ~ 6Rj and substructure at one-third this spacing 
might be due to current driven instability. Neither of the toroidally magnetized adiabatic 
jet simulations show evidence for a magnetically confined "nose cone" (Clarke, Norman, & 
Burns 1986; Lind et al. 1989). 



4.2. Magnetized Cooling Jets 

Images from three cooling jet simulations, a purely hydrodynamical model (fi£f), a 
weak field model (Pi 1 ), and a strong field model (Pu l 4 ) with magnetic field configuration 
identical to that used in the three adiabatic jet simulations described in §4.1 are shown 
in Figure 4. Figure 5 shows the toroidal magnetic field strength at the same times for 
the two magnetized jet models. Note only the final 50-R, of the computational domain is 
shown in both cases. The purely hydrodynamical jet simulation fi£f is directly comparable 
to the results presented in Stone & Norman (1993; 1994); the primary difference being the 
use here of the coronal cooling curve of Dalgarno & McCray (1972), as opposed to the 
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non-equilibrium cooling rates used previously. 

Unlike the adiabatic jets discussed above, the addition of a toroidal magnetic field to 
radiatively cooling jets leads to noticeable differences relative to the non- magnetized (3^ 
simulation. For example, in f3± l the jet beam has expanded in width, and the cocoon is 
much more uniform and smaller relative to (3^. Jet expansion in (5f l has resulted in a 50% 
decrease in the toroidal magnetic field strength in the jet beam. No nose cone is evident in 
either the or (3^ 1 simulations. 

The /3{y 4 simulation is strongly influenced by expansion of the jet resulting from the 

non-equilibrium initial configuration at the inlet. Figure 6 plots radial profiles of the density, 
axial velocity, temperature, and toroidal magnetic field strength at an axial distance of 68.5 
Rj from the jet nozzle in all three simulations. From the figure, the structure of the fi^h 
jet can be seen to consist of an expanded spine with ij pa 3.5R., surrounded by a sheath of 
thickness ~ 3.5Rj. The magnetic profile shows an approximate linear increase from the axis 
to 2Rj and plateau to 4R., with B^, pa 21 fiG. At the same time, the jet temperature drops 
from 2,000 K on axis to a low of ~ 100 K in the outer portion of the jet, and then jumps 
up to ~ 2,000 K in the cocoon. 

The stronger magnetic field along with jet expansion have resulted in significant 
morphological changes in jet and cocoon structure beyond an axial distance of about 
70 Rj in Py l A . At this distance the jet velocity suddenly decreases by more than a 

factor of two at an "upstream shock" and a nose cone of length pa 25Rj is "pushed" 
ahead of a "downstream shock". The nose cone drives a more conventional jet 
terminal and bow shock. The temporal evolution leading to this structure is shown in 
Figure 7, which plots the logarithm of the density at 5 times during the simulation, 
(0.4,0.8, 1.2, 1.6, 1.8) x t pa (61, 122, 182,243,274) yr. Initially the jet remains well 
collimated. Strong cooling at the head of the jet results in most of the shocked ambient 
and jet gas collecting in a dense shell ahead of the Mach disk. Once the jet has propagated 
about 20i?j, it begins to expand due to lack of magnetohydrostatic balance. Consequently, 
the head of the jet grows laterally, and becomes filamentary. Most of the shock processed 
gas is confined to a loosely defined "nose cone" ahead of the Mach disk. Detailed structure 
at the head of the jet is shown further in Figure 8, which plots the logarithm of the 
temperature, the divergence of the velocity (an indicator of shocks), and velocity vectors at 
the end of the simulation. Comparison of the temperature and velocity divergence plots 
clearly shows that the highest temperatures are located in filaments immediately behind 
strong shocks in the flow, such as the Mach disk and outer bow shock. Very little non-axial 
motion is observed in the velocity vectors. 



4.3. Cooling Jets with Different Magnetic Field Profiles 

Three additional magnetic field profiles have been studied for steady radiatively cooling 
jets. One simulation, fifj, has a toroidal force-free field in the jet with B<p oc 1/r. In this 
simulation the maximum field strength m = 236.3 jj,G (fl m = 0.062) at r m /Rj = 0.05 
(the innermost computational zone) and the magnetic field is set to zero at the jet surface. 
Within the force-free jet, temperature and thermal pressure are constant. In a second 
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simulation, (3^, the magnetic and pressure profiles given by equations (5) and (7) are 
used with f3 m = 0.25 (-B<^ m = 117.7 fiG) at r m /Rj = 0.2. Unlike the /3w 4 simulation, in 
this simulation a true equilibrium pressure profile can be achieved and outside the toroidal 
field maximum the field configuration is not too different from the force free simulation 
Pfj-. Finally, a third simulation, includes an axial magnetic field along with a toroidal 
magnetic field. In this simulation a constant axial magnetic field in the jet and in the 
ambient medium of strength B z = 58.8 [iG {j3 z = 1) has been added to a toroidal magnetic 
field profile given by equation (5) with, r m /Rj = 0.8 and f3 m = 0.25 (B^^ = 117.7 fiG). In 
this simulation, the pressure profile in the jet is similar to that in the (3y\ simulation and 
is not in equilibrium. Figure 9 plots the logarithm of the density at the final time for each 
of these simulations (along with the (5y A simulation shown in Figures 7 & 8), and Figure 

10 plots the corresponding toroidal magnetic field strength. Radial profiles of a variety of 
quantities in the jet beam for each simulation at an axial distance of 68.5Rj are given in 
Figure 6. 

In P^j the jet cools radiatively and the axial velocity remains constant out to ~ 90 Kj 
where the velocity abruptly drops by about 20% at a terminal "upstream" shock. A modest 
nose cone of material about 9 R, in length precedes a "downstream" shock and drives a 
terminal and bow shock. The jet expands as the initial jet pressure is above the radiatively 
cooled cocoon pressure at the inlet. The radial structure (Figure 6) is similar in dimension 
to that found in the (3^ and (3f l simulations, although the cocoon size is slightly less than 
in the (3f l simulation. 

The f3^ m simulation serves as an intermediate case between the force-free (3ff simulation 
and the non force-free f3y A simulation, and can also be compared to the (5^ 1 simulation. 

Inside the jet when r > 0.2r.,- the magnetic field profile is very similar to the force-free 
case. At the inlet the jet is overpressured relative to the radiatively cooled cocoon but jet 
expansion does not begin until an axial distance of about 40 Rj. Interior to this distance 
nine well defined pressure oscillations with a spacing of ~ 3.3 Rj are observed at jet 
radii r < 0.15r\,-, and probably are the result of a current driven pinch mode. Pressure 
fluctuations attributable to this pinch mode are on the order of ± 2% — 3% and these 
fluctuations disappear as the jet begins to expand. At an axial distance of 68.5 Rj (Figure 
6) the jet spine and sheath are reduced in size by about 25% relative to the force-free 
configuration. The jet is in approximate pressure balance with the cocoon at this point. We 
note that the jet spine contains a relatively high pressure and temperature central region 
confined by the magnetic field. The cocoon is also reduced in size by about 25% relative to 
the force-free configuration. A modest nose cone about 10 Rj long leads the jet in this case 
after a velocity drop of 15% in an "upstream" shock. 

Finally in the (5^1 simulation we consider the effect of adding an axial magnetic field 
to the toroidal magnetic field configuration used in f3y A . The axial field imposed over the 
entire computational grid has served to ameliorate the large jet expansion observed in /3jy 4 
and at an axial distance of 68.5 Rj (see Figure 6) the jet spine ij ~ 3Rj and sheath of 
thickness ~ Rj together equal the jet spine width in f3^j A . This width is similar to that seen 
for the jet and in this simulation the jet is at lower pressure than the hotter cocoon. 
Toroidal magnetic field resides within the jet spine and sheath, while the axial magnetic 
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field has been reduced within the jet by jet expansion but increased by compression in the 
outer part of the cocoon. This additional magnetic field and accompanying high sheath and 
higher cocoon density has helped to confine the jet. The presence of axial magnetic field 
has almost but not quite eliminated the extensive nose cone observed in the j3y A simulation. 



5. Pulsed Jets 

Some of the most prominent structural features in highly collimated jets associated 
with protostellar objects are the bright emission knots in the flow close to the source 
(see Reipurth 1997 for a recent review of the observations). It has been suggested that 
symmetric pinch modes of the Kelvin-Helmholtz (K-H) instability might give rise to 
such knots (Biihrke, Mundt, & Ray 1988), which has motivated detailed hydrodynamical 
modeling of the nonlinear evolution of symmetric K-H modes in cooling jets (Massaglia et 
al. 1992; Bodo et al. 1994; Rossi et al. 1997; Downes & Ray 1998). However the large proper 
motion observed in some cases (Eisloffel & Mundt 1994), and direct kinematic evidence 
provided by spectroscopic observations (Reipurth 1997) have shown that some knots are 
associated with velocity variability in the outflow. Thus, while the nonlinear stage of 
nonaxisymmetric modes of the K-H instability might still prove relevant to understanding 
the helical structures or wiggles observed in some protostellar jets (Hardee & Stone 1997; 
Stone, Xu, & Hardee 1997), most of the internal knot structure in the jets is probably 
produced by flow variability in the jet beam. 

The hydrodynamics of velocity-variable, or pulsed, jets has been widely studied in the 
literature (see Raga 1993 for a review). In an important early contribution, Raga & Kofman 
(1992) demonstrated that a sinusoidal velocity variation at the jet inlet steepens into a 
sawtooth pattern downstream. One- dimensional (Hartigan & Raymond 1993), followed 
by multidimensional hydrodynamical simulations (Stone & Norman 1993; de Gouveia Dal 
Pino & Benz 1994; Biro & Raga 1994; Biro 1996; Suttner et al. 1997) have shown that 
the velocity discontinuities in the sawtooth pattern consist of a shock pair: an upstream 
shock decelerating high velocity gas as it collides with the pulse, and a downstream shock 
sweeping up low velocity material ahead of the pulse. These shock pairs move apart at 
a constant rate causing the pulse width to grow linearly in time up to some asymptotic 
value. In the absence of magnetic fields this evolution can be understood by a simple 
analytic model (Falle & Raga 1993). Recently, there has been considerable effort in 
understanding how such "internal working surfaces" (Raga et al. 1990) affect the dynamics 
and observed properties of jets. In this section, we present the results of an investigation of 
the magnetohydrodynamics (MHD) of radiatively cooling pulsed jets. 

Parameters for the simulations discussed below are listed in Table 2. Because shock 
compression on the downstream and upstream sides of pulses generated by velocity 
fluctuation is a function of jet radius, we include in Table 2 area weighted averages as, for 
example, 

(M jtms ) = ^J q J M j:ms {r)rdr . (13) 
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5.1. Variation in Magnetic Field Strength 

We have computed a purely hydro dynamic model with (3$ = oo, hereafter denoted as 

(3^, to serve as a benchmark. Two simulations were performed using a toroidal magnetic 
field profile in the jet beam given by equation 5 in Paper II, with r m = 0.9Rj. The first 
model has a peak field strength of B<p >m = 58.8/iG corresponding to an equipartition field, 
i-e. Pfam = 8Trpo/B^ m = 1, hereafter denoted as (3[ . The second model has a peak field 
strength B^ m = 117.6/iG corresponding to a strong field where = 0.25, hereafter 
denoted as These profiles are identical to those adopted in the steady jet simulations 
(3i l , and /f/ 4 . 

In Figure 11 we show images of the logarithm of the density at the final time in the three 
simulations. The model is directly comparable to previous results of hydrodynamical 
simulations, e.g., Stone & Norman (1993). Despite the use of a different numerical method 
(a PPM algorithm was used in this earlier work, instead of the ZEUS code used here) and 
lower numerical resolution, there is excellent agreement in the overall structure of the jet. 
The steepening of the sinusoidal pulses into thin, dense sheets by z ~ 10/i,- is evident. 
Thereafter, the pulses are bounded by two shocks. Significant mass flux into the cocoon 
from the high pressure shocked gas within the pulse is evident as smooth flows. Each pulse 
is seen to drive a shock into the cocoon, so that the cocoon is much broader than steady jet 
models using the same parameters. The structure at the head of the jet is complex as pulses 
begin to merge with one another, and interact with the ambient medium. Note that the 
pulse width increases linearly with distance from the nozzle. For the purely hydrodynamic 
jet, the rate of increase in size is set by the balance between radial mass loss into the cocoon 
from the pulse "surface" and the mass flux into the pulse through isothermal shocks (Falle 
& Raga 1993). 

As the magnetic field strength is increased (models /3f and /J^ 4 shown in the lower two 

panels of Figure 11), the gross properties of the pulsed jet remain the same. Two features of 
the pulses are clearly modified by the magnetic field. Firstly, the postshock density within 
the pulses on the axis increases with the field strength. Secondly, the rate of increase in 
the pulse width increases with increasing magnetic field strength. While the presence of 
hoop stresses associated with the toroidal field clearly influences the amount of material 
ejected from the pulses into the cocoon, the maximum radial distance that material is 
ejected into the cocoon does not decrease with field strength. In all three simulations, the 
wings generated by the pulses extend to a radial distance of ~ 5Rj. Figure 12 plots the 
logarithm of the toroidal magnetic field strength in (3[ and f3y A . The field acts as a tracer 

of jet material since the ambient gas is unmagnetized. Compression of the magnetic field 
by the internal shocks associated with the pulses, and ejection of material and the magnetic 
field into the cocoon is evident. Note this ejection occurs over a larger area in /3^ 4 as might 
be expected given the larger pulse width. 

In order to illustrate the variation of properties of the pulses with magnetic field 
strength, we plot in Figure 13 the density, pressure, axial velocity, and toroidal magnetic 
field strength along an axial slice through the center of all three jets at time t = 1.2r - 
note that this time is not the same as that used in Figures 11 and 12. Only a limited 
segment of the axial domain between 20Rj and 45Rj is shown for clarity. The plot shows 
three pulses developing within this region, one centered near 23Rj, the next at 31 Rj, and 
the last at 38Rj. In both cases the magnetic field strength declines outwards from the 
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origin over the region shown in Figure 13 by a factor 2-3. Several trends are clearly 
evident in the plot. Note the increase in width of the pulses with axial distance. The 
strongly magnetized model shows the widest pulses. Interestingly, the gas pressure is nearly 
identical within all three jets. The axial density is much higher in the strongly magnetized 
pulsed jet, increasing from 5 times (in the first pulse) to 20 times (in the last pulse) larger 
than the hydrodynamic pulsed jet. The axial velocity shows the classic sawtooth pattern 
expected for a sinusoidally pulsed jet (Raga & Kofman 1992), modified by the shock pair 
which appears at each step (e.g., Stone & Norman 1993). The toroidal field strength follows 
the density profile closely as would be expected in ideal MHD when the magnetic field is 
parallel to the shock interface. 

In all three jets, the pulse width increases linearly with axial distance. At z ~ 50-R-,-, 
the pulse widths are about 4Rj, 5Rj, and 7Rj in (3^, (3[, and (3[, 4 , respectively. To zeroth 
order the linear growth rate in the pulse width scales with 1/ (M Jims ) rather than the sound 
speed (Falle & Raga 1993) (see Table 2). In (3[ and (3^ the pulses remain separated along 
the entire length of the jet, so that the upstream and downstream shock pair is evident 
nearly to the jet head. On the other hand, the spreading of the pulse width is so large in 
the f3y 4 model that the pulses merge before reaching the head of the jet, leading to the 
disappearance of the upstream shock - leaving only a sequence of downstream shocks - at 
z > 60Rj. 

The observed differences in pulse density, radial density profile and pulse width are the 
result of hoop stresses which confine some of the shocked jet material within the pulses 
to near the axis, and the result of the lack of radial pressure equilibrium in the material 
after it passes through the upstream and downstream shocks into the pulse. Recall that 
the unshocked material in (3[ is in radial pressure equilibrium. In an isothermal shock the 
toroidal magnetic pressure increases oc p 2 whereas thermal pressure increases oc p. Thus, in 
general, the toroidal plasma (3 behind the shocks should decrease relative to the toroidal 
plasma (3 in the unshocked jet material by up to oc p _1 . Of course, the toroidal plasma (3 in 
the unshocked and shocked pulse material is a function of radius and the relative importance 
of thermal versus magnetic pressure is further modified by material and magnetic ejection 
from the pulse. However, in j3^ 4 the area weighted plasma (3 in a typical pulse is ((3^) ^ 2 

(compared to (fy) ~ 8 in the unshocked jet material), and in this case magnetic pressure 
can play a significant role in supporting the pulses against ram pressure. 

Considerable toroidal field accumulates at the head of the jet into a nose-cone-like 
structure; this is most evident in the very strong field case. For our parameters typical 
magnetic field strengths in the pulses at z > 20Rj are about 50pG and 75/zG in (3[ and 
f3y 4 , respectively, and the field strength declines as the distance from the inlet increases. 
Between the pulses the magnetic field strength is reduced from these values by over a factor 
of 10. Near the head of the jet in the complex interactions between shocked jet material 
ejected from the pulses and shocked ambient gas, field strengths are as high as 190pG in 
(3[ , and 570pG in /3f /4 . 
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5.2. Radial Variation in Magnetic Field and Pulse Structure 

We have studied the effect of varying the geometry of the magnetic field in the jet 
beam using three different simulations. The first, is the model discussed above in 
which -B^m = 58.8/iG at r = 0.9Rj. The second, /3^ m , is another equipartition toroidal 
magnetic field model but with B^m = 58.8/iG at r = 0.2Rj. The third, (3% z , is a model with 
B<f> tm = 58.8/^G at r = 0.8Rj, and which contains a uniform axial field of B z = 58.8/iG, 
with f3 z = 8irp /B z = 1 added throughout the computational domain. In all three of these 
simulations, the jet beam is initially in exact radial equilibrium. 

Images of the logarithm of the density at the final time in each of these simulations are 
shown in Figure 14. A comparison between f3f and f3^ m shows that when the magnetic field 
is peaked close to the axis, the pulses widen non-uniformly, growing in width more rapidly 
near the axis of the jet than at the surface. The size and shape of the wings of the pulses is 
similar in both models. A strong axial magnetic field (model f3^ z ) helps to confine material 
in the pulses and reduces the amount of material ejected to the cocoon. 

Images of the logarithm of the toroidal magnetic field strength in (3^ m , and (3% z 

are shown in Figure 15. In , considerable toroidal magnetic field has been injected into 
the cocoon. Since the field tracks jet material, this indicates a large amount of shocked 
jet gas is ejected from the pulses. On the other hand, in f3^ m the strongest toroidal field 
is found close to the axis, and there is much less magnetic field in the cocoon, indicating 
reduced ejection of shocked jet material from the pulses. In (3% z the toroidal magnetic field 
is more confined towards the axis indicating that considerably less material is ejected from 
the pulses into the cocoon, and the cocoon is smaller in this simulation than in the other 
two simulations. 

In Figure 16 we plot various quantities as a function of radius for pulses at an axial 
distance z/Rj = 48, 53, 51, 47.5 and 51 at the final time in simulations (3[, f3y A , (3^ m 

and /3% z , respectively. These locations are the centers of last pulse which is non-interacting 
with its neighbors at the final time in each simulation. Note that (3[/ 4 and (3^ m result in 
similar axially peaked density pulse profiles (although the axial density in /5^ 4 is four times 

that in /3^ m ). Other magnetic field configurations are less axially peaked but are still very 
different from the "top hat" profile of the hydrodynamical simulations (solid line in each 
panel). In all cases the magnetized pulses are significantly colder on the axis, typically by 
more than a factor of 2, than the hydrodynamic pulse. The hydrodynamic pulse has a "top 
hat" temperature profile with T vu \ se ps 1, 300 K. Expansion velocities evaluated at r = 0.8-Rj 
lie within ±15% of ctj (dj = a /3.16) and are relatively independent of magnetic field 
strength or configuration. Thus, we might expect an equilibrium pulse width (mass into 
the pulse through upstream and downstream shocks balanced by transverse mass ejection) 
to scale approximately proportional to the density at the r = 0.8Rj surface, i.e., with the 
mass ejection rate per unit surface area evaluated at r = 0.8Rj. 

At least approximately we would expect the pulse expansion rate and ultimate width 
in /3i and j3^ z to be comparable but twice that (and in f3[, 4 and /3^ m comparable but about 

four times that) observed in the hydrodynamic case. Our simulations do in fact show an 
increase in the rate at which pulses widen consistent with this analysis, although the rate 
of this increase is not as fast as this simple argument would predict. 
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From Figure 16, the maximum radial ejection velocity from the pulses occurs at a radial 
distance much larger than the pulse surface at r = 0.8Rj. It is affected by the magnetic 
field configuration and strength, with the highest velocities observed for the strongest 
toroidal magnetic field case, fl[u, and lowest velocities observed when an axial magnetic 

field is present, f3^ z . Changes to the unshocked toroidal plasma (3 in the pulses is also 
dependent on the magnetic field configuration and strength, and the radial variation in 
shown in Figure 16 is much different from that at the inlet (see Table 2). If we characterize 
the effects of the magnetic field by the area weighted toroidal plasma (3 then we find that 
{(Impulse ~ (fit) jet in cases (3[ and @% z , whereas ((3<j,)p Ulse « ((3<f,)j e t in cases (3y A and (3^ m . 
Note that this properly groups the magnetic simulations with their observed spreading 
rates. However, only in case [3[/ 4 where (P^pui se ~ 1.5 can there be significant magnetic 

pressure effects. We conclude that the differences in pulse width spreading rate are the 
result of the radial density profile in the pulse induced by radial pressure equilibrium with 
13$, and clearly the density profile is significantly modified even when the area weighted 
toroidal plasma (5 is much greater than one. 



6. Summary 



6.1. Steady Jets 

Although inclusion of a weak magnetic field (weak in the sense that the flow is highly 
supermagnetosonic) has little effect on the propagation of adiabatic, overdense jets (e.g. 
Figures 2 and 3), such fields have strong effects on cooling, overdense jets (e.g. Figures 4, 5, 
9 and 10). The difference is caused by the much larger compression of jet and ambient gas 
that occurs in cooling shocks. This compression amplifies the magnetic field until in some 
cases the magnetic pressure is comparable to the gas pressure in shocked gas, a situation 
which never occurs in the adiabatic jet simulations. 

Our simulations reveal at least three effects of the magnetic field on steady, cooling 
jets. Firstly, the field affects the structure and fragmentation of dense sheets of cooled gas 
(e.g. Figures 4 and 9). This is in part because magnetic pressure can support cooling gas 
against compression, and so limit the density and width of cooling shells. Alternatively, 
magnetic stresses can inhibit the radial flow of gas away from the axis, leading to higher 
densities in some cases. Both of these effects are evident in comparison of models (3^, 

and /3f/ 4 (Figure 4). The structure of the magnetic field is important to its overall 
effect, for example, helical magnetic fields appears to produce smoother fluctuations in 
the cocoon than purely toroidal fields of the same strength (see Figure 9). Much of this 
may be attributed to the difference the field has on the radial structure of the jet, either 
through expansion or compression of the jet beam by pressure and hoop stresses. Jets 
which undergo radial expansion have a much larger cross section near their tip, and drive 
larger, more fragmented shells into the ambient gas. Jets which are confined or compressed 
by the magnetic field tend to have narrower and smoother cocoons. 

Secondly, the presence of a magnetic field can, in some cases, result in the production of 
a dense, plug of shocked jet material ahead of the Mach disk (a "nose cone"). The formation 



-18- 



of nose cones has been noted previously in simulations of magnetized extragalactic jets 
(Clarke et al. 1986). However, we find the formation of nose cones in cooling jets is sensitive 
to the geometry of the field: only purely toroidal fields which peak near the surface of the 
jet form nose cones, simply because this configuration maximizes the hoop stresses which 
confine the nose cone. Moreover, nose cones are likely to unstable in 3D which will limit 
their relevance to real protostellar jets. 

Finally, the presence of magnetic fields affects the stability of the jet beam. We do not 
observe any evidence for the MHD K-H instability in our simulations, nor do we expect 
to given that the linear growth rates for such modes are small for the supermagnetosonic 
jets studied here. Interestingly, we observe axial pressure fluctuations in a few simulations 
where the jet beam is unstable to current driven pinch modes (Begelman 1998). However, 
our simulations reveal such modes saturate at low amplitude in the nonlinear regime, so 
that it is unlikely they will have much relevance to the internal knots observed in most 
protostellar jets. 



6.2. Pulsed Jets 

Our simulations reveal that weak magnetic fields have a number of effects on 
time-variable (pulsed) cooling jets. Such pulsing results in the formation of dense knots of 
shocked jet material in the jet beam. The inclusion of a toroidal magnetic field affects the 
rate at which such knots grow by inhibiting the radial flow of material out of the knots and 
into the cocoon. This results in much higher densities in the pulses as the field strength is 
increased (Figure 13). The increase of the density in the cooling gas in the shocked pulses 
as the magnetic field strength is increased contradicts the expectations of simple planar 
shock theory: this result emphasizes the importance of multidimensional effects in MHD 
flows. 

A second effect of a non-uniform magnetic field is the introduction of radial structure 
into the density and pressure of the pulses, even if the jet is uniform at the nozzle. Jets with 
strong toroidal fields have radial density profiles in the pulses which are strongly peaked 
toward the axis. This profile is introduced by radial variations in the magnetic pressure and 
hoop stresses in the shocked pulses, which confines material near the axis. It is unlikely 
that real protostellar jets have a uniform ("top-hat") profile: it has been shown that the 
structure of the cocoon can be affected by the profile of the density and velocity in the jet 
beam through in hydrodynamical simulations (Suttner et al 1997). Our results reinforce 
the importance of understanding the radial profile of the outflow to interpret the structures 
observed. 

Finally, even though toroidal fields affect the internal radial structure of the pulses, we 
do not observe much effect on the structure of the cocoon due to the magnetic field. Radial 
profiles of various quantities in the pulses (Figure 16) show that the maximum axial velocity 
of material ejected from the pulses into the cocoon is relatively insensitive to the field 
strength. Thus, the size of the cocoon is not significantly decreased in the magnetized jets. 
The radial ejection velocity is small (a few times the sound speed in the unperturbed jet 
beam), and it would require a delicate force balance to completely eliminate radial outflow. 
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7. Conclusions 

We have studied the propagation of both steady and time- variable ("pulsed") 
protostellar jets using numerical MHD simulations. Although we focus our attention on 
models in which the jet beam contains a purely toroidal magnetic field peaked near its 
surface, we have also studied the effect of varying both the strength and geometry of the 
field on the dynamics. 

We find that even a weak magnetic field (B < 60/xG) in the jet beam can lead to 
important effects on the structure and dynamics of steady jets (FRJN-C; CGH; CG). In 
particular, such fields can alter the density structure and fragmentation of dense shells 
formed in cooling jets. However, the details of the effects depend sensitively on the geometry 
of the field. For example, while magnetic pressure can limit compression in such shells, 
resulting in lower densities and less fragmented shells, radial hoop stresses associated with 
purely toroidal fields can confine shocked gas towards the axis, resulting in higher densities 
there. In some cases, hoop stresses lead to the formation of nose cones ahead of the 
Mach disk. We see no evidence for significant structure induced by the Kelvin-Helmholtz 
instability (such modes were seen for the parameters and field geometry adopted by CG). 
We do see evidence for current driven pinch modes but the induced pressure pulses in the 
jet interior do not appear strong enough to explain the emission knots that are observed 
in protostellar jets (see also CG). Thus, models invoking jet velocity fluctuation appear to 
remain the most viable explanation for the knots observed in protostellar outflows. We do 
not see any evidence for disruption of the jet by pinch modes. 

In the case of pulsed cooling jets, we find the primary effect of a toroidal magnetic 
field is to confine shocked jet material to the axis, preventing it from begin ejected into the 
cocoon, and leading to higher postshock densities in the pulses in comparison to purely 
hydrodynamic models. This result is in contrast to the expectation of planar radiative 
shock models, in which the addition of a magnetic field leads to lower postshock densities. 
Our results indicate it is important to account for multidimensional effects in the study of 
magnetized cooling jets. Toroidal confinement also leads to radial variation of quantities in 
the pulses, even if the density and velocity initially are constant with radius. The radial 
variation of, e.g. the density, in the pulses is large enough (a factor of 80 in a strongly 
magnetized jet) that it will likely affect the computation of the resulting emission properties. 

While there are many uncertainties in the magnetic field strength and topology implied 
by the observations, our present results make it clear that MHD models of protostellar 
jets need to be seriously investigated. This is particularly true since we find that average 
values of the magnetic field which correspond to a plasma (5 much larger than one (and 
which therefore one would infer to be too weak to be important) can still make a significant 
difference in the dynamics and physical conditions associated with the jet, and therefore its 
emission properties as well. It is thus important to understand the asymptotic structure of 
the field produced by the mechanism which drives the outflow. 
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8. Figure Captions 

Fig. 1— Radial profile of the toroidal magnetic field (solid line) and thermal pressure 
(dashed line) in the jet beam. 

Fig. 2. — Images of the logarithm of the density at the final time in adiabatic overdense 
jets with no magnetic field - - (top panel), and weak - j3f d - (middle panel) or strong 
-fiijl - (bottom panel) toroidal field. A linear color scale between logn = 0.97 and 4.32 is 
used, where n is the particle number density in cm -3 . The maximum density in each case is 
logn = 4.05, 4.32 and 3.48, respectively. The initial density of the jet and ambient medium 
are logn, = 3 and logn a = 2. The axes are labelled in units of the initial jet radius; note 
only the last 50Rj of the computation is shown in each case. For typical protostellar jet 
parameters, r = 152 yrs and lOOii, = 0.08 pc. 

Fig. 3. — Images of the logarithm of for the magnetized jets shown in Figure 2. A 
linear color map between \og(B ( p/ ^/4tv) = —2 and 1.1 is used. The maximum values are 
\og(B ( p/\^ATr) = 0.54 and 1.1, respectively. 

Fig. 4. — Images of the logarithm of the density at the final time in cooling overdense jets 
with no magnetic field - (3^ - (top panel), and weak - (3± l - (middle panel) or strong - /3^ l A 
- (bottom panel) toroidal field. A linear color scale between logn = 0.30 and 5.6 is used. 
The maximum density in each case is logn = 4.2, 5.6 and 5.4, respectively. 

Fig. 5. — Images of the logarithm of B^ for the magnetized jets shown in Figure 4. A 
linear color map between \og{B^/ yf^n) = —2 and 1.4 is used. The maximum values are 
\og{B ( j > /\f^K) = 1.16 and 1.4, respectively. 

Fig. 6. — Radial profiles of density, axial velocity, temperature, and toroidal magnetic field 
in cooling steady jet simulations at z ~ Q8.5Rj from (3™ (solid), {3f l (dotted), /3^L (long 
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dash), p^j (dot long dash), (3^ m (dot short dash), and (3^ l z (short dash). The density, axial 
velocity, and temperature are scaled relative to the ambient density, p a , the fiducial speed, 
ao, and temperature T = 10, 000 K. 

Fig. 7.— Time history of the logarithm of the density during the propagation of an 
overdense, cooling, strongly magnetized jet. A linear color scale between logn = —0.06 
and 5.4 is used. 

Fig. 8. — Logarithm of the temperature (top panel), velocity divergence (middle panel), and 
velocity vectors (bottom panel) near the head of the overdense, cooling, strongly magnetized 
fliL jet, at time 1.8r ~ 274 yrs. A linear color map is used for the temperature between 
logT = 1.87 and 5.8. Only negative velocity divergence values are shown, sharp filaments 
indicate the location of shocks in the flow. The maximum length of velocity vectors is scaled 
to 350 km s -1 . 

Fig. 9. — Images of the logarithm of the density at the final time in cooling overdense jets. 
Top panel: a strong toroidal field - (3y 4 - which peaks at r = 0.9Rj. Second panel: a strong 
toroidal field - (3^ - which peaks at r = 0.2Rj. Third panel: a strong toroidal field - (3^ l z 

- which peaks at r = 0.9Rj combined with a weak axial field. Bottom panel: a force-free 
toroidal field - Pfj. A linear color scale between logn = 0.36 and 5.4 is used. The maximum 
density in each case is logn = 5.36, 4.93, 4.33 and 5.36, respectively. From top to bottom 
the jets are shown at times of 1.8r , 1.8r , 1.7r , and 1.7t respectively 

Fig. 10. — Images of the logarithm of for the magnetized jets shown in Figure 9. A 
linear color map between \og{B^/ \f^a) = —2 and 1.5 is used. The maximum values are 
log(S^/V47r) = 1.37, 1.22, 1.22 and 1.49, respectively. 

Fig. 11. — Images of the logarithm of the density at the final time for the pulsed jets: @£> 

- no magnetic field - (top panel), (3[ - equipartition toroidal field - (middle panel), and 
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(3y 4 - strong toroidal field - (bottom panel). A linear color scale between logn = — 1.4 and 
5.3 is used, where n is the number density in cm -3 . The maximum density in each case is 
logn = 4.46, 4.65, and 5.29 respectively. The initial density of the jet and ambient medium 
is logn, = 3 and logn a = 2 respectively. The axes are labelled in units of the initial jet 
radius; note the vertical scale is expanded by a factor of two for clarity. 

Fig. 12. — Images of the logarithm of for the magnetized pulsed jets (3[ and (3y A shown 
in Figure 11. A linear color map between log^^/V^r) = —2 and 1.3 is used. The maximum 
values are \og(B ( p/ 1 \J\~k) = 1.2 and 1.3, respectively. 

Fig. 13. — Axial profiles along the axis (r = 0) of density, pressure, velocity, and magnetic 
field for pulsed jets (solid line), f3[ (dotted line) and /3^ 4 (dashed line) at t = 1.2r . 
Only the region between 20Rj and 45i?^ from the jet nozzle is shown. 

Fig. 14. — Images of the logarithm of the density at the final time for the pulsed jets: f3^ A 
- strong toroidal field peaked at r m = 0.9Rj - (top panel), - equipartition toroidal 
field peaked at r m = 0.2Rj - (middle panel) and (3^ z - equipartition toroidal field peaked 
at r m = 0.8i?j and a uniform axial field - (bottom panel). A linear color scale between 
logn = —2.1 and 4.9 is used. The maximum density in each case is logn = 4.65, 4.86, and 
4.35 respectively. Note the vertical scale is expanded by a factor of two for clarity. From top 
to bottom the jets are shown at times of 1.6, 1.6, and 1.7r respectively. 

Fig. 15. — Images of the logarithm of for the magnetized pulsed jets /3f/ 4 , (3^ m , and /3£ z 
shown in Figure 14. A linear color map between log(-E>^/ '\/A~k) = —2 and 1.3 is used. The 
maximum values are log(-B^/v^47r) = 1.2, 1.3, and 1.1, respectively. 

Fig. 16. — Radial profiles of density, radial velocity, thermal pressure, and toroidal plasma 
P in pulses at z ~ 50Rj from (3^ (solid line), (3[ (dotted line), /3^ m (dashed-dotted line), f3^ z 
(short dashed line), and (3[< A (long dashed line). The density, radial velocity, and thermal 
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pressure are scaled relative to the ambient density, p a , the fiducial speed a , and the inlet 
thermal pressure at the jet edge, p . 
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TABLE 1 

Simulation Parameters For Steady Jets 



Simulation 


M [a) 

J,IILt> 


M {b) 

J ,// lb 




P, {a) /Po 


f vn 1 -^j 






1 oo 


63.2 


63.2 


20 


1 








i l 


43.9 


48.3 


20 


2 


0.9 


1 




oAd 


99 4 


9^ ^ 


90 


c 

O 


n q 








63.2 


63.2 


283 


1 








p? 1 


43.9 


48.3 


283 


2 


0.9 


1 




oCl 
Pl/A 


22.4 


25.5 


283 


8 


0.9 


0.25 




oCl 
Pff 


14.1 


30.1 


283 


1 


0.05 


0.062 




PRm 


19.2 


42.0 


283 


11 


0.2 


0.25 




qCl ( d ) 

Pd>z 


20.9 


24.5 


18.3 


8 


0.8 


0.25 


1 



a on jet axis in first computational zone 
b linear average 

c 3 = 1 coresponds to = 58.8^(7 
d 10 zones/Rj 



-27- 



TABLE 2 

Simulation Parameters For Pulsed Jets 



Simulation 


M (a) 

lv± j,ms 










I'm/ Rj 


3 (c) 






/£ 


63.2 


63.2 


63.2 


283 


1 








oo 




43.9 


50.5 


68.9 


283 


2 


0.9 


1 




6.7 




22.4 


27.3 


46.4 


283 


8 


0.9 


0.25 




8.2 




33.6 


57.2 


60.7 


283 


3.5 


0.2 


1 




16.4 


fa 


34.5 


38.3 


64.6 


18.3 


2.2 


0.8 


1 


1 


10.2 



a on jet axis in first computational zone 

b area weighted average 

c f3 = 1 coresponds to B — 58.8/iG 
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